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Abstract 

Cooperation is a very common, yet not fully-understood phenomenon in natural and human 
systems. The introduction of a network within the population is known to affect the outcome of 
cooperative dynamics, allowing for the survival of cooperation in adverse scenarios. Recently, the 
introduction of multiplex networks has yet again modified the expectations for the outcome of the 
Prisoner’s Dilemma game, compared to the monoplex case. However, much remains unstudied 
regarding other social dilemmas on multiplex, as well as the unexplored microscopic underpinnings 
of it. In this paper, we systematically study the evolution of cooperation in all four games in the 
T — S plane on multiplex. More importantly, we find some remarkable and previously unknown 
features in the microscopic organization of the strategies, that are responsible for the important 
differences between cooperative dynamics in monoplex and multiplex. Specifically, we find that in 
the stationary state, there are individuals that play the same strategy in all layers (coherent), and 
others that don’t (incoherent). This second group of players is responsible for the surprising fact 
of a non full-cooperation in the Harmony Game on multiplex, never observed before, as well as a 
higher-than-expected cooperation rates in some regions of the other three social dilemmas. 
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1. INTRODUCTION 


Cooperation is a ubiquitous and yet not fully-understood phenomenon in Nature: from 
humans that cooperate to build complex societies to animals like wolves that hunt in packs 
in order to catch preys larger than they are, or meerkats that watch out for predators in turn 
while the rest of the colony feeds. Even small microorganism cooperate to survive in hostile 
environments. For instance, the Dictyostelium discoideumu, usually a solitary amoeba, when 
starves it associates with others in order to form a multicellular slug for the sake of survival. 
Explaining how cooperation has emerged and has resisted against more selhsh behaviours 
is one of the biggest challenges in natural and social sciences. From a mathematical point 
of view, the problem of cooperation within a population can be studied using Evolutionary 
Game Theory mEi. There are multiple mechanisms proposed to explain the evolution of 
cooperation, such as reputation, kin selection, network reciprocity or punishment |1^E]. On 
the other hand, outstanding experimental efforts have been made in the last few years [Ml] 
to try to understand how actual humans behave when confronted with social dilemmas in a 
formal Game Theory environment. 

We focus here on the impact of the structure of the network of interactions among individ¬ 
uals on the outcomes of the cooperation dynamics. The study of networks, their properties 
and dynamics, has experimented a huge advance in the last few decades, empowered by the 
technological advances that enable the acquisition of real data about interactions between 
individuals from social networks [laiisi, mobile communication networks na or collabo¬ 
rations between scientihc authors [I5|. There is a vast literature on the evolution of coop¬ 
eration on complex networks [T61 - 1T8] . studying aspects ranging from the effect of network 
topology on cooperation [19] to network growth driven by cooperation dynamics [201 1^ - 
and other spatial and temporal effects [22] that offer insights on how cooperation can evolve 
and survive in different scenarios. 

An innovative way of representing multiple types of social interactions in one single struc¬ 
ture is the use of multiplex networks [23H2S], see Fig. which have been already successfully 
applied to the study of disease spreading IZH and diffusion dynamics [2H] (for a complete 
review look at H). Multiplex networks are interesting in this held, because many social 
interactions can be understood as a combination of interactions at different, independent 
levels, each one representing a different social scenario such as family, friends, coworkers. 
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FIG. 1. Example of a multiplex network with 3 layers, 5 nodes per layer and 5 links in each layer. 
The color of the nodes represents the strategy played in that layer, red for cooperators, blue for 
defectors. Their size is proportional to their global payoff. 

etc. An individual’s behaviour can be different in each level, but it is ultimately conditioned 
by all of them. Some work has been done to understand the evolution of the Prisoner’s 
Dilemma game on multiplex networks |30], exploring different coupled evolutionary games 
using a interdependent networks [31]. The impact of the degree correlations among lay¬ 
ers [32] on the outcome of social dilemmas have also been studied on 2-layer network, where 
one layer was used for the accumulation of payoffs and the other for strategy updating. 
There are also works that explore the problem of cooperation on coupled networks [33], 
and even optimizing the interdependence between them via coevolution [M]. However, the 
evolution of cooperation on top of multiplex networks with any number of layers hasn’t been 
systematically studied for all four social dilemmas. 

The objective of this paper is, on the one hand, to provide an exhaustive analysis of the 
resilience and propagation of cooperation in the main four social dilemmas in Game Theory 
literature, studying the average levels of cooperation, payoff distribution, and dependence 
on the initial fraction of cooperation, as a function of the number of layers of the multi¬ 
plex. More importantly, we will focus on analyzing the previously unexplored microscopic 
behaviour of individuals across layers. 

This work is organized as follows. In Section we dehne the model we have used in 
this work. Section [^ contains our hndings on the density of cooperators for each one of the 
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proposed scenarios. Then we turn our attention to the microscopic behaviour of individuals 
across different layers. Finally, a summary and conclusions can be found in Section]^ 


2. MODEL 

We will focus on two-strategy social dilemmas. If we assume that each player in the 
system can either cooperate (C) or defect (D), a game can be dehned according to its payoff 
matrix: 

C D 

C IR S 
d\t P 

Where R represents the reward obtained by a cooperator playing against another cooper¬ 
ator, S is the sucker payoff obtained by a cooperator when it plays against a defector, the 
temptation payoff, T, is the payoff received by a defector when his opponent is a cooperator, 
and hnally, P represents the payoff obtained by a defectors which engages with another 
defector. 

Traditionally the values of R and P are hxed to i? = 1 and P = 0 in order to provide 
a hxed scale for the game payoffs [35l [36]. Applying this constraint, it turns out that 
the selection of the remaining parameters T and S enables the dehnition of several games 
according to their evolutionary stability. Thus, if P > S' > P and R > T > P the game 
is the harmony game 133. The hnal state of a population playing this game will be total 
cooperation, regardless of the initial fraction of cooperators. Prisoner’s dilemma [SHHlQ], 
T > R > P > S, represents the opposite situation, and the population evolves towards 
total defection regardless of the initial conditions (although all players would be better off 
cooperating, hence the dilemma). A classical example of a coordination game, the stag-hunt 
game 1111112], is represented when the payoff values respect the order R > T > P > S, the 
output of this game will be either total defection or total cooperation, depending on the 
initial conditions. Finally, an anti-coordination game, the Hawk-Dove mm, takes place if 
the payoff values follows T > R > S > P, where the hnal state will be a population made 
of both cooperators and defectors. 

The players sit on the nodes of a multiplex network of L layers. Each node is present 
in all layers, but in general, they have different connectivity in each layer. Every layer, A, 
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in the multiplex network is a connected and homogeneous Erdds-Renyi (ER) network, with 
the same number of edges E and nodes iV, and equal degree distribution, and the multiplex 
network is generated avoiding degree correlations between layers. Each layer is represented 
by an adjacency matrix where = 1 if nodes i and j are connected in that layer, and 
Rfj = 0 otherwise. That representation enables the dehnition of the degree of node i in layer 
I as kl = its global degree in the multiplex as Ki = X]r=i 

Each round of the game is divided in two phases: payoff recollection and strategy up¬ 
date. Each node i can choose to play one of the two strategies, cooperation or defection, 
independently in each layer of the network and at every time step, sf(f). Within a specihc 
payoff matrix, the node i’s strategy determines the payoff, pi, that it obtains in a layer I 
when it plays against all its kl neighbors. The total payoff of node i can be easily calculated 
as Pi = J2f=oPi- round, each player can change the strategy in one of its 

layers, sf, using the Replicator-like Rule: A node chooses a layer of the multiplex, £r, with 
uniform probability. Then it chooses with uniform probability one of its kl’’ neighbors, >, 
in that layer. If P, < Pj^ and sf 7 ^ the probability that node i changes its strategy in 
layer is given by: 


_ P3rit)-P^{t) _ 

max(Arj, Kj^) ■ (max(l, T) — min(0, S)) 

It is important to notice that the update rule uses global information about the players: 
global degree and global payoff (that is, added up over all layers), in order to update the 
strategy of any particular layer. That is the way our model shares information between 
layers and relies in the social nature of layers’ interdependency |30]: each player only has 
information about the strategy of its neighbour in their same layer (but not in those layers 
where they are not connected). However, it knows its neighbor’s total benehts, and it makes 
the simplifying assumption that it is using the same strategy in every layer. As we will see 
later on, this fact has a profound impact on the outcomes of the dynamics, compared to the 
monoplex scenario. 

At the end of each time step the density of cooperators can be computed for each layer 
and for the entire multiplex using: 



N 
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e=i i=i 


(3) 
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3. RESULTS 


To ascertain the outcome of the cooperative dynamics for the different games on multiplex 
networks, we will start by studying the stationary level of cooperation in the system, then 
we will study the effect of the initial fraction of cooperators, and hnally, we will move to 
analyzing in detail the microscopic organization of cooperation for individuals across different 
layers. 

The results are obtained for a range of values of T G [0,2] and S G [—1,1] that dehnes 
the T — S plane. The simulation runs on a multiplex network that has N = 1000 nodes and 
E = 3000 edges per layer distributed according an Erdos-Renyi degree distribution with (k). 
For each possible pair of values of the game parameters the simulation runs 1 x 10® time 
steps, that is the transient time to needed by the algorithm to generally reach a stationary 
state (we further discuss the matter of convergence time in the Supplementary Information, 


Section SI). After this time the algorithm runs for another = 2 x 10^ time steps. All 
the quantities of interest are averaged over this second period of time. The experiments are 
repeated and averaged over / = 64 different networks and initializations in order to gain 
statistical conhdence. The initial fraction of cooperators, cq, is distributed randomly in each 
layer. We focus here on the case cq = 0.5, although we have also explored other values (see 


Supplementary Information, Section S5). 

Density of cooperators. The stationary average value of cooperation is dehned ac¬ 
cording to the following: 


^ i=l t=to 

In Fig. 1^ we present the average stationary value of cooperation when cq = 0.5. We 
observe that our results for the monoplex case (left) are consistent with those obtained by 
Roca et al. [22] for this kind of networks. The results for multiplex networks show a large 
increase of the areas where both strategies coexist (that is, the areas in the plane that 
separate total cooperation from total defection). However, this coexistence is of a different 
nature depending on the evolutionary stability of the particular game (or quadrant), as we 
explain below. 

The Stag Hunt game has an unstable evolutionary equilibrium with mixed populations. 
This means that, when there is a structure, the population will evolve towards total co- 
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Monoplex (L = 1) Multiplex (L = 5) Multiplex (L = 10) 

1.000 0.734 0.932 0.696 0.910 0.661 



T 


FIG. 2. Asymptotic density of cooperators (c) for networks with different number of layers (L = 1 
on the left, L = 5 in the middle, L = 10 on the right). The plane T — 5 is divided into four major 
regions that correspond to the four games under study: the upper-left area is the Harmony Game, 
the upper-right is the Snow Drift, Stag-Hunt is in the lower-left, and the Prisoner’s Dilemma is 
in the lower-right. The average asymptotic density of cooperators for each one of the games is 
also indicated, as a numerical value, next to the corresponding quadrant. See Supplementary 


Information, Section S5, for the corresponding results for other values of the initial fractions of 
cooperation. 


operation or total defection depending on the initial population and type of structure of 
the network (due to this fact, the standard deviation of the (c) is large in that transition 


area, see Supplementary Information, Section SI, for details). For the monoplex we have a 
very narrow transition area between total cooperation or defection populations (left panel 
in Fig. 1^. This transition region widens with the number of layers, enabling the coexistent 
of both strategies in a larger portion of the game parameter space. The explanation of such 
behaviour can be found in the inter-layer dynamics: it is more likely that a cooperator or a 
defector node resists in hostile environments in a particular layer, because its htness is not 
evaluated in just that layer, but also in the other layers where, due to its strategy or its 
topological conhguration, the node might have better performance. The Stag Hunt game, 
where the maximum payoff possible is obtained when a cooperator plays against another 
cooperator, favors specially the resilience of cooperators nodes when the temptation value 
is low: a cooperator node i in layer that has a big payoff Pi has higher probability of 
spreading its strategy to its defector neighbours in thus increasing its payoff. This in- 
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crease will propagate to the other layers, making the strategies of the player more robust 
against invasion. Playing defection in layer ir when temptation value is small, does not have 
a big effect in the global payoff of the node. As a consequence, in this particular game the 
multiplex structure increases specially the resilience of cooperators, thus the average density 
of cooperators in this game quadrant shows an statistically signihcant increase as we keep 
adding layers to the structure (Mann-Whitnet U test, a = 0.05). 

In the Prisoner’s dilemma game, defection dominates cooperation. Related papers [22] 
show that for ER networks using Replicator rule when temptation and sucker payoffs are 
not too large, cooperation can survive forming groups of cooperative clusters, thus resisting 
against the initial attempt of invasion by defectors, and then spread through the population. 
Our results for the monoplex are consistent with those. For the multiplex, we observe how 
the transition region between all-cooperator and all-defector situations is larger than for the 
monoplex, as in the case of Stag Hunt game. It is worth noticing that regions where we have 
all-cooperator populations in the monoplex, are not necessarily all-cooperator areas in its 
multiplex counterpart. This happens because the formation of cooperative clusters in one 
layer will also increase the htness of these nodes in the other layers regardless of the strategy 
used in these other layers. And this can lead to a reinforcement of defector strategies due 
to the inter-layer dynamics, increasing their survival rate. This inter-layer dynamics will 
led to a widening of the transition area that enables survival of cooperators in areas where 
they are not present in the monoplex scenario. If we take into account the whole Prisoner’s 
Dilemma quadrant, the conclusions are the same that in the Stag Hunt game: a statistically 
signihcant increase in the average density of cooperators occurs as we increase the number 
of layers. 

The Snow Drift game has a stable equilibrium in mixed populations: it is an anti¬ 
coordination game. Previous works show that for ER networks there are some regions 
in the plane T — S for which this game converges to single-strategy populations. For lower 
values of the temptation these regions are prone to cooperation. In multiplex networks 
however, single strategy regions are less common and mixed populations are the rule. That 
happens by the same inter-layer dynamics that we have explained earlier: the impact of a 
cooperator’s benehts on the other layers of the multiplex structure. This entails a signihcant 
reduction on the average fraction of cooperators from 0.734 in the monoplex to 0.661 in the 
10-layer multiplex for this quadrant. 


Finally, the Harmony game has cooperation as its dominant strategy. For single-layer 
ER networks with Replicator update rule, Roca et ah [22] reported that the whole quadrant 
ends up in an all-cooperator conhguration. However, in the case of multiplex scenarios, 
the average fraction of cooperators decreases significantly as we keep adding layers to the 
system: 0.932 for L = 5 and 0.910 for L = 10. This increasing resilience of defection can 
be explained as a consequence of the multiplex topology and the lack of degree correlations 
between layers: due to the payoff accumulated by an individual acting as cooperator in some 
layers, defector nodes can resist against cooperators in other layers. 

We can mathematically prove that defectors can survive and be stable in the Harmony 
game on ER multiplex networks by analyzing the simplest situation: let’s assume a multiplex 
structure with L layers. In one single layer (for simplicity we assume it will be the hrst one) 
we have one single node playing as defector, but it plays as cooperator in all the other L — 1 
layers. There are no more defectors anywhere in the system. This node’s connectivity in 
layer a is ka, and, recalling that R = 1 and P = 0, the total payoff of that node that is 
defecting in one single layer is given by: 


L 

Pd = Tkl ka (5) 

a=2 

The payoff of any of the node’s neighbors (note that all of them play as cooperators), 
with a degree k'^ in layer a, is: 

L L 

a=«-l) + S + ^A.; = ^A^; + S-l (6) 

a=2 j=l 

Thus, in order to survive as a defector in layer a, the following inequality must be fulhlled 
for each of the node’s neighbours: 


Pd > Pc ( 7 ) 

L L 

Tfci + ^ ^ fc' + ^ - 1 (8) 

i=2 j=l 

We can estimate both a soft and a hard limit for the previous inequality. As a soft 
limit, and assuming we have independent, uncorrelated Erdos-Renyi layers in our multiplex 
network, we can approximate every k^ by {k) and get: 
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(T + L- l)(fc) >S-1 + L{k) 

(9) 

S -1 

(10) 

{k)>T_, 


On the one hand, a hard limit for the condition can be calculated by approximating ki 
by fcmax for the cooperator neighbours: 


{T + L-l){k) >S-l + Lk^ 
\ . S' — 1 + Lk„ 

(k) > 


T + L 


( 11 ) 

( 12 ) 


On the other hand, we can calculate the probability of this topological situation hap¬ 
pening. First of all we have to dehne what is the probability of a node i to have degree fc, 
P{X = k). In our model, and in order to avoid the non-negligible effect of unconnected 
nodes, we impose a minimum connectivity, fcmin- To get a more accurate approximation of 
our degree distribution we take into account this minimum: 


. {X = k) = -—- 

As it has been stated previously, the payoff of cooperators against cooperators is propor¬ 
tional to their degree, since we set R = 1: in this example we use L — 1 full cooperative 
layers, so the payoff obtained in this layers is proportional to the degree distribution of the 
aggregate network of this L — 1 layers. Moreover, the payoff distribution of the nodes that 
play cooperation in all layers is proportional to the aggregation of all layers, L. Imposing 
that we do not have inter-layer degree correlation, the degree distribution of the aggregated 
networks can be modeled using the convolution of the single layer degree distributions. 



Pl 


rsj 



(14) 


The probability that a topological conhguration that enables the fulhlment of the payoff 
conditions specihed by Eq. ([^ exists, is given by: 

OO OO 

Psurvivai= Pk^,^{X = k,) Pl_,{X = q) ■ Pl{X < [q + k, ■ T - S + (15) 

^l=^inin 
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Multiplex (L = 5) Multiplex (L = 10) Multiplex (L = 100) 



0.0 0.5 1.0 0.0 0.5 1.0 0.0 0.5 1.0 

T 


FIG. 3. Probability of a defector surviving in the Harmony Game for 5 layers (left), 10 layers 


(middle) and 100 layers (right), calculated according to Eq. (15). The individual layers are ER 
with (k) = 3. 


where q is the payoff obtained by the defector node playing as a cooperator in L — 1 layers. 
With that information, an upper bound for the aggregated degree of the defector’s neigh¬ 
bours can be dehned as [q + ki ■ T — S + 1\, and if all the neighbours have an aggregated 
degree below this upper bound, the defector can survive. It is worth noticing that the upper 
bound for the degree of a cooperator is a discretization of payoff values that involve S and 
T. This means that the survival probability of a defector only changes when the relation 
between S and T changes by an amount large enough. 

The expression for the degree distribution probability function is for an Erdos-Renyi net¬ 
work, assuming that we have a restriction for the minimum degree, so the degree distribution 
follows a Poisson distribution given by: 


k 


n„n(V = k) = 


A*'e 


(16) 


«■ (1 - #) 

In Fig. we show the probability of a defector surviving in a full-cooperative population. 


calculated numerically using Eq. (15). We observe that this probability increases naturally 
with T, because this is the payoff that a defector obtains against a cooperator, but it is 
only slightly dependent of the payoff of a cooperator against a defector, S. The number 
of layers has a huge impact on this probability: as the number of layers increases, the 
probability becomes more uniform in the T — S plane, increasing in general. This can 
be explained by the relative contribution to the accumulated payoff that comes from layer 
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1 (the layer where the defector survives): the more layers are added to the system, the 
smaller this relative contribution. For a large number of layers, this implies that the values 
of S and T (that determine the payoff) are less important in the probability of a defector 


persisting in the system. For networks with a higher mean degree (see Fig. S4 in the 
Supplementary Information), however, the chances of a defector surviving are lower: if the 
number of neighbours of the defector node is higher, then the probability that one of them 
has more payoff than him is also higher, thus the defector will tend to imitate the neighbour’s 
behaviour (or in other words, his chances of survival will decrease). 

Coherent Cooperation. Prompted by the topological conhgurations described earlier, 
we can now dehne a “coherent cooperator” as a node that, at a given instant of time, plays 
as cooperator in all L layers of the system. Similarly, we can dehne a “coherent defector” 
as a node that, at a given instant of time, plays as defector in all L layers of the system. 
Finally, those individuals that are neither coherent cooperators nor coherent defectors will be 
called “incoherent” individuals. This new terms introduced here should not be mistaken for 
the concepts “pure cooperators”, “pure defectors” and “huctuating individuals” introduced 
in [19], which implied a temporal consistency of the agents’ strategies. Also, we want to 
stress that a incoherent individual as dehned here, is clearly different from the concept of 
a mixed population, that refers simply to a set of both strategies, coexisting together in a 
population. Moreover, we have to take into account that a coherent behaviour is not trivial 
nor easily reachable, due to the fact that our simulations start with all mixed populations 
(randomly distributed and uncorrelated strategies in all layers), so the dynamics that leads 
to coherence is specially interesting to study. 

In Fig. I^we show the fraction of coherent cooperators (left column), coherent defectors 
(middle column) and incoherent individuals (right column) for 5 layers (top row) and 10 
layers (bottom row). The formation of coherent cooperators is particularly complicated, 
and it is interesting to notice that even in the Harmony game there is a low fraction of 
them (except for a small area around the extreme case of T = 0 and S' = 1). In the other 
quadrants, the fraction is very small (in particular, the Prisoner’s Dilemma presents basically 
no coherent cooperation). This implies that most of the cooperation shown by the system 
comes from incoherent individuals. We also observe that the fraction of coherent cooperators 
decreases quickly with the number of layers for any game. As we have said, the origin of such 
results resides in the fact that a defector takes advantage of its own cooperative behaviour 
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Coherent Cooperators Coherent Defectors Incoherent Players 


0.697 0.259 0.000 0.033 0.303 0.709 



T 


1.00 


0.75 


■ 0.50 


i 0.25 

I 


0.00 


FIG. 4. Average density of coherent cooperators (left column), coherent defectors (middle column) 
and incoherent individuals (right column) for networks with 5 layers (top row) and 10 layers (bottom 
row). The average density of the corresponding type of individuals is also provided for each one of 
the quadrants (upper-left is the Harmony Game, upper-right is the Snow Drift, Stag-Hunt is the 
lower-left, and the Prisoner’s Dilemma in the lower-right). 

in other layers, specially in regions of the T — S plane prone to cooperation. 

Conversely, regarding the fraction of coherent defectors, we observe that their presence is 
very strong in most of the Prisoner’s Dilemma region and part of the Stag-Hunt area, and 
they decrease only slightly when increasing the number of layers from 5 to 10. This fact is 
easy to understand: the resilience of a cooperator in a hostile environment is based basically 
in how he performs as cooperator, the advantage of playing as defector in other layers is 
practically zero because in a large defector population the contribution to the payoff of a 
defector that plays against a defector is zero, P = 0. Thus, in these regions, the survival 
rate of cooperation does not improve by playing as defector strategy in other layers. 

Regarding incoherent individuals, we observe that they are very prevalent for all games 
(except for the extreme area of Harmony around (T = 0, S' = 1), where cooperation is very 
prohtable, and the bottom-half area of the hard Prisoner’s Dilemma where cooperation is 
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extremely expensive). Incoherent individuals contribute significantly to the average density 
of cooperation in a large central area of the T—S plane, particularly in the areas that separate 


full-cooperation from full-defection (See also Fig. S5 in the Supplementary Information for 
a detailed description of the fraction of incoherent individuals playing as cooperators). This 
area of prevalent incoherent individuals increases with the number of layers or, in other 
words, it gets harder and harder to be a coherent strategist as the number of layers increases. 

Fig. 1^ also confirms what we showed analytically earlier: defection can survive in the 
Harmony game, as long as the individual defecting in a particular layer is a incoherent 
individual; it plays as cooperator in other layers and obtains enough payoff from them to 


avoid having to switch strategies (see also Fig. S5 in Supplementary Information for further 
detail on the payoff of cooperators and defectors). 

Interestingly enough, in Fig. we can observe how coherent players of opposite types 
do not coexist in the same population. Another important point is where coherent players 
can coexist with incoherent players. The area where coherent cooperators interact with 
incoherent players is wide and gets wider as we keep adding layers to the multiplex. However, 
the area of coexistence of coherent defectors and incoherent players is very narrowed and is 
only slightly affected when layers are added to the structure. This means that the coherent 
defection is a very dominant strategy that almost forbids the existence of any other kind of 
players. 


4. DISCUSSION 

In this paper we have presented a systematic and comprehensive analysis of the outcomes 
of cooperation dynamics on ER multiplex networks for the four games on the T — S plane, 
when using the Replicator updating rule, comparing our results with those already known 
for the case of the games on monoplex. Also, we have analyzed the microscopic behavior of 
the nodes, and coined the terms of coherent cooperator, coherent defector and incoherent 
player. 

In particular, we have found that the stationary distribution of cooperation in the plane 
T — S becomes less sharp as more layers are added. In the monoplex case there is a very 
narrow area that separates all-cooperator from all-defector areas for the Stag Hunt and 
Prisoner’s Dilemma games, but in the multiplex scenario we find that it becomes a wider 
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region, with intermediate values of cooperation. We also find that the region of all-defectors 
shrinks as the number of layers increases. As a counter-effect though, we hnd a slight 
decrease in the value of cooperation (even in the quadrant of the Harmony game), from total 
cooperation to values around 90%. These results are consistent with and generalize those 
found by PI: the introduction of a multiplex structure in the population helps promote 
cooperation in regions of the parameter space in which it can not survive in the monoplex 
scenario, at the expense of a moderate decrease of cooperation in those where traditionally it 
was very high. We explored the microscopic underpinnings for these phenomena, previously 
observed but unexplained in the aforementioned paper. 

Thus, regarding the microscopic behavior of the nodes, we have found that in general and 
at a given time step, there are three types of individuals: those coherently acting as cooper¬ 
ators in all layers, those acting as coherent defectors, and a group of incoherent individuals, 
that play as cooperators in some layers and as defectors in others. The existence of this third 
incoherent group is at the root of the explanation of the survival of defection in the Harmony 
Game for a multilayered network, and it is also responsible for a large part of the cooperation 
in the central areas of the T — S plane, where cooperation is lower in a monoplex. Also, we 
have analyzed how this three types of players interact among them, concluding that there 
are plenty of interaction between incoherent and coherent cooperators, fewer interactions 
between incoherent and coherent cooperators, and practically no interaction between both 
types of coherent players. Moreover, this is a very plausible social scenario: some people 
may behave consistently in all their types of interactions (for example at work, at home, with 
friends,etc) either cooperating or defecting, and some other may choose different strategies 
for different layers (for example, cooperate with family and defect at work). We have found 
that an the fraction of incoherent players increases with the number of layers increases, 
which means that as the number of contexts where the a players interact increases, it gets 
harder to maintain a coherence behaviour in all of them. Regarding the dependence with 
the initial fraction of cooperation, we found that our system behaves consistently with what 
was found for the monoplex network, and the effect of adding more layers is preserved or 
even increased with increasing initial fraction of cooperators. 

To summarize, the introduction of multiplex networks not only is a more realistic repre¬ 
sentation of social systems, allowing for more sophisticated individual behaviours, but as it 
has been shown in other context too, it has a profound effect on the dynamics developing 
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on top of them. 
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Supplementary Information: Strategical incoherence regulates 
cooperation in social dilemmas on multiplex networks 


SI. CONVERGENCE 


We study the system’s convergence to the stationary state. It is well known that, in 
general, the time evolution of cooperation on a monoplex network, for a value of the pa¬ 
rameters that allows the survival of at least some cooperation, usually follows a curve that 
initially decreases moderately, while cooperation rearranges itself from the random initial 
conditions into a more favorable setting (either in one or multiple clusters), and then there 
is a new increase, followed by the achievement of the stationary state. In general, this whole 
transient time is relatively short (typically of 1 — 2 x 10^) for a size of 1 — 4 x 10^ nodes. 
However, it hasn’t been explored in detail until now the convergence process for the four 
games in the T — S plane on multiplex networks. 

In order to evaluate such convergence, we £t the last time steps of the evolution to 
a linear trend, c{t) = a + pt using the QR decomposition method. Then we use the slope 
of the fitted model to compute the variation of the density of cooperators every 1000 time 
steps, S = 1000 ■ (3. Thus, a near-zero value of this metrics indicates that the system has 
reached the stationary state, while a positive value would indicate that the average level of 
cooperation is still increasing in the system at that time, and vice versa. Figure [ST] shows 
how every point of the T — S plane performs on our measure of convergence during the last 
time steps of the simulation. 

Monoplex networks (left column in Figure [ST| seem to reach the stationary state according 
to our convergence criteria for every point of the plane T — S and independently of the initial 
fraction of cooperators: the slope is in general smaller than 10“"^. We observe, however, a 
small amount of stochastic noise for some regions of Snow-Drift and Prisoner’s Dilemma 
games, where our measure indicates that the stationary is not fully reached. Nonetheless, 


we will show in Section that this noise is just an effect of the large fluctuations in the 
number of cooperators when the stationary is reached. 

In multiplex networks, on the other hand, there is a non-negligible area where convergence 
is not reached (red areas in the central and right panels in the hrst row in Figure [ST|. In 
the most extreme cases, where the slope of the linear model (3 is largest, our measurements 
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FIG. SI. Convergence to the stationary state measured as the average variation in the fraction of 
cooperators during 1000 time steps, measured using the slope of a linear model fitted at the end of 
steps of the simulation. The numbers in each quadrant represent the mean convergence value 
(the slope of the fit) for each one of the four games (upper-left is the Harmony Game, upper-right 
is the Snow Drift, Stag-Hunt is the lower-left, and the Prisoner’s Dilemma in the lower-right). 
In the different rows we show the information for several values of initial fraction of cooperators 
(co = 0.25, Co = 0.50, Co = 0.75), while the different columns correspond to 1, 5 and 10 layers, 
respectively. 


indicates an increment of the cooperators of about a 0.1% every 1000 time steps. That could 
seem a smaill increase in the fraction of cooperators, however if the evolutionary process 
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FIG. S2. Example of the time evolution of cooperation for the point in the plane T — S with 
maximum fluctuations values (S' = 0 and T = 1.35). The shadowed area for each plot represents 
two standard deviations over the residuals of the I iterations at each time step. It will be used to 
compare the size of the fluctuations between monoplex and multiplex networks. The grey vertical 
area corresponds to the interval + ^ 7 ] where the measures shown in all panel figures in this 

paper are computed. 


were to run for a very long period of time, the increase could be significant. 

To better illustrate this difference in the path to stability for monoplex vs. multiplex 


networks, we show in Figure S2 the time evolution of the level of cooperation, (c) for a single 
simulation progresses (monoplex plotted in red and two multiplex networks with different 
number of layers represented in green and blue), for one point in the plane T — S. This 
particular point has been picked as an extreme case, for having the maximum fluctuation 
values in the entire T — S plane (see Figurej^and Sectionj^for further detail). We clearly 
observe that, while the time required for the monoplex system to achieve the stationary 
state is around 1 — 2 x 10^, for the multiplex networks it can be at least one order of 
magnitude larger, and it increases with the number of layers, too. However, it is important 
to remember that this example shown here is a very extreme case, while the convergence 
process in multiplex is in general faster for regions of the plane that are far away from the 
transition area (or areas where the hnal state is close to an all-cooperation or all-defection). 

To understand the reason for such an increase in the convergence time for multiplex with 
respect to monoplex (at least in some regions of the T — S plane), one has to pay attention 
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to which areas are more reluctant to reach stability. Such regions correspond again to the 
transition areas between those that end up in total cooperation and those that end up in 
total defection. In the Stag Hunt quadrant, the game has an unstable equilibrium with 
mixed population, which means that the game will tend to converge to total cooperation or 
total defection as happens in the monoplex network. However the multiplex structure of the 
network changes that outcome, as we described in the Results Section. In these structures, 
the transition region is larger than in the monoplex, and is in this transition region where 
the convergence is hard. The analysis of how the fraction of cooperators has an effect on the 
convergence gives us an insight about what is happening. We have already stated that the 
interlayer dynamics has an important role in the survival rate of defectors and cooperators. 
If we look at the multiplex columns of the Figure [STj we can observe how the convergence is 
strongly affected by the initial fraction of cooperators. On the one hand, if the initial fraction 
of cooperators is small, they will need more time to reach the equilibrium because they have 
to hght against a larger fraction of defectors that benehts from the interlayer dynamics. On 
the other hand, a larger initial number of cooperators implies that the defectors will need 
more time to reach an stable conhguration. However, the presence of a large number of 
initial cooperators has less impact on the convergence; which is easily understood, given the 
fact that defectors get more proht from cooperating in other layers than the opposite case. 

Similar conclusions could be reached for the other games, taking into account that the 
transition regions between full-cooperation and full-defection are different in nature, for 
instance in the Snow-Drift this region is wider. Thus, we can see the effect of non-convergence 
is diluted across the Snow-Drift quadrant. 


S2. ANALYSIS OF FLUCTUATIONS 

We turn our attention now to the fluctuations of the system in the stationary state. In the 
case of these four games on a monoplex network, it is well known that the level of cooperation 
in the stationary state fluctuates around a well-dehned average value due to the effect of 
both the topological structure of the network and the nature of the Replicator updating 
rule. We propose a measure in order to quantify these fluctuations and later compare them 
with the cases of multiplex networks. For each one of the I repetitions of the experiment 
we £t a linear model to the final t-y time steps of the simulation, Ci{t) = ai + I3it. We also 
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need to take into account a possible non-zero slope in the measure of fluctuations (see in 


Section SI), so we average the Mean Square Error between the data from the simulations 


and the predictions of the linear model for the I iterations, calculated as: 


^ / h+to 

( = jY, = T— E E ('=•(*) -'5(*))' 


i=l 


(SI) 


^ i=l t=to 


The results for monoplex and multiplex networks are displayed in Figure S3 For the 
monoplex case, the simulations show small fluctuations in the quadrants of Stag-Hunt and 
Harmony Game. However for the Snow-drift on monoplex, the results display a zone where 
fluctuations are larger than in the rest of the plane T — S. That can be attributed to 
the nature of the game: it has an evolutionary stable equilibrium with mixed populations, 
so a consensus where both strategies coexist has to be reached. To achieve this objective 
some nodes have to alternate their strategies. These changes, due the topological features 
of the network, can lead to a cascade effect of changes in a large portion of the network; 
the equilibrium gets disturbed, and a new equilibrium has to be reached again. This causes 
the relatively large fluctuations that we measure. It is worth noticing that in the area 
of mutual coexistence of strategies the fluctuations are larger where the temptation and 
sucker payoffs are not far from the payoffs of mutual cooperation and mutual defection. It 
is also noteworthy that, even when the Prisoner’s Dilemma quadrant presents very small 
fluctuations in general, it does show a small but very signihcant spot near the line of weak 
Prisoner’s Dilemma, where they are large. Again, this corresponds to the area of competition 
between Cooperation and Defection, where each of the strategies accounts roughly for half 
the population. 

The introduction of multiplex networks has an enormous effect on the fluctuations. The 
fluctuations are again in the region of coexistence of strategies, however, in the case of 5- 
layer multiplex the fluctuations are much smaller than in the monoplex case. The results 
in the 10-layer multiplex display an even larger reduction in the measure of fluctuations 
(compare also the three example curves shown in Figure [S^. The nature of such reduction 
from monoplex to multiplex is to be found in the interlayer dynamics. Each layer is pushed 
to reach an stable equilibrium where both strategies can coexist, nonetheless the shared 
information between the layers establishes a way to constrict the range of the fluctuations. 
The change of strategy of a node in one layers is not conditioned by its performance in that 
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FIG. S3. Fluctuation of the fraction of cooperator around the fitted trend at the final time steps 


of the simulation computed as described in the equation (SI), for each pair of values S and T. The 
number in each quadrant represents the average value for each one of the four games (upper-left is 
the Harmony Game, upper-right is the Snow Drift, Stag-Hunt is the lower-left, and the Prisoner’s 
Dilemma in the lower-right). The results are provided by 3 different initial conditions 


single layer, but by its global performance in the entire multiplex structure. That makes 
the system more robust to fluctuating nodes, at the expense of convergence time to the 
stationary equilibrium. We observe that, both for monoplex and multiplex structures, the 
initial fraction of cooperators, cq, barely has any observable influence on the size of the 
fluctuations. 
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FIG. S4. Probability of a defector surviving in the Harmony Game for 5 layers (left), 10 layers 
(middle) and 100 layers (right) for an average connectivity of (k) = 20, calculated according to 
equation 15 from the main text. 


It is worth mentioning that the fluctuations shown in Figure are calculated with a 


modihed version of equation (SI), as follows: we have to divide the time range in different 
slices in the interest of realizing local accurate measures of the fluctuations in each slice, so 
we £x the size of the time window to = 1000. For each time slice we fit a linear model to 
each of the I runs of our simulation, then we compute the residuals as the difference between 
the linear model and the data from the simulations. We have I residuals that measure the 
size of the fluctuations at each time step, and we plot a range corresponding to twice their 
standard deviation to provide information about the size of the fluctuations at each time 
step. 


S3. PROBABILITY OF A DEFECTOR SURVIVING IN THE HARMONY GAME 
FOR HIGHER AVERAGE DEGREE 

In Figure we show the probability of a defector surviving in a full-cooperative pop¬ 
ulation, calculated numerically using equation 15, for the case of a higher average degree, 
(k) = 20, than in the main text. The main impact of an increased value of average degree 
is a signihcant decrease of the probability, for any number of layers, or values of S and T. 
The effects discussed in the main text remain for this case too, but attenuated (note that 
the range of values for the probability are smaller in this case). In general, the probability 
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increases with T, but it is only slightly dependent of S. As the number of layers increases, 
the probability becomes more uniform in the S — T plane, increasing in general. 

S4. PERCENTAGE OF COOPERATION AMONG MIXED INDIVIDUALS, AND 
PAYOFF OF COOPERATORS AND DEFECTORS 

In this hnal section, we address the analysis of the percentage of cooperation among mixed 
individuals, as well as the payoff obtained by both cooperators and defectors. Regarding the 
former, we observe that the percentage of mixed individuals playing as cooperators is very 
high in the Harmony game, and in the upper diagonal of the Stag-Hunt game, as well as 
the upper diagonal of the Snow-Drift. In the Prisoner’s Dilemma game, however, it is zero 
except for a small region near the weak limit, when cooperation is relatively inexpensive. 
This general situation gets emphasize by the increasing of the number of layers. 

Regarding the latter, we observe that the average payoff among cooperators in obviously 
the highest in the Harmony game, and upper diagonals of both Stag-Hung and Snow-Drift, 
and it is zero in the hardest, bottom diagonal of the Prisoner’s Dilemma game, with a wide 
transition area of intermediate values separating both regions. Moreover, this description 
seems to be independent of the number of layers in the system. Finally, the only regions 
where defection gets a moderate payoff are within the Snow Drift game, while it is near zero 
anywhere else. This picture is also independent of the number of layers. 


S5. INITIAL FRACTION OF COOPERATORS. 


For completeness, in Fig. we show the stationary average fraction of cooperation for 
the four-game plane and various numbers of layers, for three different initial fractions of 
cooperators (upper row is cq = 0.25, middle row is cq = 0.5 and bottom row is cq = 0.75). 
We will briefly discuss now the differences between the previously explained case of Cq = 0.5, 
and the other two scenarios. 

We observe that, for a given game quadrant and a given number of layers, increasing the 
initial fraction of cooperation has in general a positive but moderate impact on the stationary 
fraction of cooperation, specially in the Stag Hunt and Prisoner’s Dilemma games. In the 
former one, we have an unstable evolutionary equilibrium in mixed populations, so the 




FIG. S5. Percentage of cooperation in mixed individuals (left column), average direct payoff 
obtained playing as cooperator (middle column) and average direct payoff obtained playing as 
defector (left column) for 5 layers (top row) and 10 layers (bottom row) multiplex, in the four 
games. The corresponding averaged values over the quadrants are also provided (upper-left is the 
Harmony Game, upper-right is the Snow Drift, Stag-Hunt is the lower-left, and the Prisoner’s 
Dilemma in the lower-right). 

change of cq has a signihcant impact on the hnal outcome. In the case of Prisoner’s dilemma 
game, an increase in the initial fraction of cooperators means an increase in the probability 
that clusters of cooperators forms. 

The effect discussed in this paper when adding layers to the system still holds or is even 
emphasized by an increased initial fraction of cooperators: the overall stationary value of 
cooperation increases with the number of layers in the Prisoner’s Dilemma Game and Stag- 
Hunt, the region of coexistence between both strategies widens for the Snow-Drift Game, 
and the Harmony game presents a small decrease of cooperation. 
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FIG. S6. Asymptotic density of cooperators (c) for networks with different number of layers (L = 1 
in the left column, L = 5 in the central column, L = 10 in the right column), and different initial 
fraction of cooperation (cq = 0.25 in the top row, cq = 0.5 in the central row, cq = 0.75 in the 
bottom row). The plane T — S' is divided into four major regions that correspond to the four games 
under study; the upper-left area is the Harmony Game, the upper-right is the Snow Drift, Stag- 
Hunt is in the lower-left, and the Prisoner’s Dilemma is in the lower-right. The average asymptotic 
density of cooperators for each one of the games is also indicated, as a numerical value, next to 
the corresponding quadrant. 
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